import numpy as np
import matplotlib.pyplot as plt

accel_A
accel_B
abfourier


#P12 - 3
#Calcule os coeficientes 𝑎𝑛 e 𝑏𝑛 de Fourier do movimento
#Faça o gráfico de 𝑎𝑛 e 𝑏𝑛 em função de 𝜔𝑛.


# Funções de aceleração
def accel_A(xA, xB):
    return -(k / m) * (xA - xA_eq) - (k_ / m) * (xA - xB)

def accel_B(xA, xB):
    return -(k / m) * (xB - xB_eq) - (k_ / m) * (xB - xA)


# Integração numérica (Euler-Cromer)
for i in range(n_steps - 1):
    aA = accel_A(xA[i], xB[i])
    aB = accel_B(xA[i], xB[i])
    
    vA[i+1] = vA[i] + aA * dt
    xA[i+1] = xA[i] + vA[i+1] * dt
    
    vB[i+1] = vB[i] + aB * dt
    xB[i+1] = xB[i] + vB[i+1] * dt



# Cálculo dos coeficientes de Fourier para xA
it0 = 0
it1 = len(t) - 1
N = 30
omega = 2 * np.pi / T

a_vals = []
b_vals = []
omega_n = []

for n in range(1, N + 1):
    af, bf = abfourier(t, xA, it0, it1, n)
    a_vals.append(af)
    b_vals.append(bf)
    omega_n.append(n * omega)





def abfourier(tp,xp,it0,it1,nf):
#
# cálculo dos coeficientes de Fourier a_nf e b_nf
#       a_nf = 2/T integral ( xp cos( nf w) ) dt   entre tp(it0) e tp(it1)
#       b_nf = 2/T integral ( xp sin( nf w) ) dt   entre tp(it0) e tp(it1)    
# integracao numerica pela aproximação trapezoidal
# input: matrizes tempo tp   (abcissas)
#                 posição xp (ordenadas) 
#       indices inicial it0
#               final   it1  (ao fim de um período)   
#       nf índice de Fourier
# output: af_bf e bf_nf  
# 
    dt=tp[1]-tp[0]
    per=tp[it1]-tp[it0]
    ome=2*np.pi/per

    s1=xp[it0]*np.cos(nf*ome*tp[it0])
    s2=xp[it1]*np.cos(nf*ome*tp[it1])
    st=xp[it0+1:it1]*np.cos(nf*ome*tp[it0+1:it1])
    soma=np.sum(st)
    
    q1=xp[it0]*np.sin(nf*ome*tp[it0])
    q2=xp[it1]*np.sin(nf*ome*tp[it1])
    qt=xp[it0+1:it1]*np.sin(nf*ome*tp[it0+1:it1])
    somq=np.sum(qt)
    
    intega=((s1+s2)/2+soma)*dt
    af=2/per*intega
    integq=((q1+q2)/2+somq)*dt
    bf=2/per*integq
    return af,bf